home *** CD-ROM | disk | FTP | other *** search
/ Info-Mac 4 / Info_Mac IV CD-ROM (Pacific HiTech Inc.)(August 1994).iso / Science / RLaB / examples / nn-test.r < prev    next >
Text File  |  1994-04-25  |  3KB  |  157 lines
  1. //
  2. // @(#) mail.r
  3. //
  4. // ([{
  5. //
  6.  
  7. // necessarily global parameters:
  8. //   fircoeff
  9. //   nfir
  10. //   nsaved
  11. //   sigsave
  12. //   firout
  13. //   desired
  14. //   iter
  15.  
  16. // initialize model filter
  17. initfir = function () {
  18.  
  19.     // model filter parameters
  20. //    fircoeff = [0.2;-0.4;1;-0.6;0.4];    // model filter coefficients
  21.     fircoeff = [1];
  22.     nfir = size(fircoeff)[1];
  23.     nsaved = 20;            // signal delay line length, must
  24.                     // be >= nfir
  25.     rand("default");
  26.     sigsave = rand(nsaved, 1);        // initialize
  27.  
  28. }
  29.  
  30. // model filter time step
  31. stepfir = function () {
  32.  
  33.     // excitation signal generation
  34.     rand("default");
  35. //    signal = (rand() > 0.5) * 0.8 + 0.1; // exp
  36. //    signal = (rand() > 0.5) * 1.8 - 0.9; // tanh
  37.     signal = (rand() > 0.5) * 1.8 - 0.9; // simple
  38.  
  39.     // filter update
  40.     sigsave = [signal; sigsave[1:(nsaved - 1)]];
  41.     firout = sigsave[1:nfir]' * fircoeff;
  42.  
  43.     // nonlinearity
  44.     // firout = tanh(3 * firout);
  45.     // firout = 2 * sin(3 * firout) + firout;
  46.  
  47.     // noise
  48.     rand("normal", 0, 0.1);
  49.     firout = firout + rand();
  50.  
  51.     // impulse noise
  52.     rand("default");
  53.     if (rand() > 0.7) {
  54.         firout = firout + 5;
  55.     }
  56.  
  57.     // what we would like to adapt to
  58.     // desired = firout;
  59.     desired = sigsave[3;];
  60.  
  61. }
  62.  
  63. bp2init = function () {
  64.  
  65.     // backprop parameters
  66.     ninput = 10;            // input vector length
  67.     input = zeros(ninput, 1);    // initialize
  68.     ndf = 0;            // number of feedback taps
  69.     nhidden = 10;
  70.     rand("default");
  71.     weight1 = (2 * rand(nhidden, ninput ) - 1) * 0.1;
  72.     theta1  = (2 * rand(nhidden, 1      ) - 1) * 0.1;
  73.     weight2 = (2 * rand(1      , nhidden) - 1) * 0.1;
  74.     theta2  = (2 * rand(1      , 1      ) - 1) * 0.1;
  75.     mu = 1.0;
  76.  
  77.     // miscellanea
  78.     niter = 20000;
  79.     skiplimit = 0;
  80.     nwsave = niter - skiplimit;
  81.     wsave = zeros(nwsave, 5);
  82.     esquare = 0;
  83.     output = 0;        // so that we can shift the first value into input[]
  84.  
  85.     initfir();
  86.  
  87. }
  88.  
  89. bp2steps = function (fromiter, toiter, skiplimit) {
  90.  
  91.     for (iter in fromiter:toiter) {
  92.  
  93.         mu = 0.01 * (niter - iter + 1) / niter;
  94.  
  95.         stepfir();
  96.  
  97.         // reconstrunction filter data generation
  98.         input = [firout; input[1:(ninput - 1)]];
  99.  
  100.         // forward pass
  101.         s1 = weight1 * input + theta1;
  102.         hidden = s1 ./ (1 + abs(s1)); // simple
  103.         s2 = weight2 * hidden + theta2;
  104.         output = s2 ./ (1 + abs(s2)); // simple
  105.  
  106.         // error vectors
  107.         e2 = (mu * (1 ./ ((1 + abs(s2)) .* (1 + abs(s2))) + 0.05)) .* (desired - output); // simple
  108.         e1 = (1 ./ ((1 + abs(s1)) .* (1 + abs(s1))) + 0.05) .* (weight2' * e2); // simple
  109.  
  110.         // weight update
  111.         weight1 = weight1 + e1 * input';
  112.         weight2 = weight2 + e2 * hidden';
  113.  
  114.         // threshold update
  115.         theta1 = theta1 + e1;
  116.         theta2 = theta2 + e2;
  117.  
  118.         err = desired - output;
  119.         esquare = 0.995 * esquare + 0.005 * (err' * err);
  120.         if (iter > skiplimit) {
  121.             wsave[(iter - skiplimit);] = [desired, firout, output, err, esquare];
  122.         }
  123.  
  124.     }
  125.  
  126. }
  127.  
  128. bp2 = function () {
  129.  
  130.     iter = 1;
  131.     toiter = 0;
  132.     while (iter < niter) {
  133.         fromiter = toiter + 1;
  134.         toiter = fromiter + 999;
  135.         if (toiter > niter) {
  136.             toiter = niter;
  137.         }
  138.  
  139.         bp2steps(fromiter, toiter, skiplimit);
  140.  
  141.         printf("Iteration %d: esquare = %g,\n", iter, esquare);
  142.         printf("   weight1^2 = %g, theta1^2 = %g,\n", trace(weight1 * weight1'), theta1' * theta1);
  143.         printf("   weight2^2 = %g, theta2^2 = %g\n", weight2 * weight2', theta2' * theta2);
  144.     }
  145.  
  146.     plotdata = [((skiplimit + 1):niter)', wsave];
  147.  
  148.         plgrid();
  149.         ptitle ( "RLaB Neural Net Example (contributed)" );
  150.         xlabel ( "" );
  151.         ylabel ( "" );
  152.     plot(plotdata[nwsave-100:nwsave;]);
  153.  
  154. }
  155.  
  156. // }])
  157.